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We have computed line emission cooling rates for the main cooling species in models of 
interstellar molecular clouds. The models are based on numerical simulations of super-sonic 
magneto-hydrodynamic (MHD) turbulence. Non-LTE radiative transfer calculations have been 
performed to properly account for the complex density and velocity structures in the MHD 
simulations. 

Three models are used. Two of the models are based on MHD simulations with different mag- 
netic field strength (one model is super-Alfvenic, while the other has equipartition of magnetic 
and kinetic energy). The third model includes the computation of self-gravity (in the super- 
Alfvenic regime of turbulence) . The density and velocity fields in the simulations are determined 
self-consistently by the dynamics of super-sonic turbulence. The models are intended to rep- 
resent molecular clouds with linear size L w 6pc and mean density (n) w 300 cm~ 3 , with the 
density exceeding 10 4 cm~ 3 in the densest cores. 

We present 12 CO, 13 CO, C 18 0, O2, 01, CI and H 2 cooling rates in isothermal clouds with 
kinetic temperatures 10-80 K. Analytical approximations are derived for the cooling rates. 

The inhomogeneity of the models reduces photon trapping and enhances the cooling in the 
densest parts of the clouds. Compared with earlier models the cooling rates are less affected by 
optical depth effects. The main effects come, however, from the density variation since cooling 
efficiency increases with density. This is very important for the cooling of the clouds as a whole 
since most cooling is provided by gas with density above the average. 

Subject headings: ISM: clouds - radio lines: ISM - ISM: molecules - ISM: structure - Radiative Transfer 
1. Introduction peratures of molecular clouds are determined by 



ABSTRACT 



Line emission by molecules and atomic species 
is the most important cooling process in interstel- 
lar clouds (Goldsmith & Langer 1978). The tem- 
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the balance between radiative cooling and various 
heating mechanisms like cosmic ray heating, for- 
mation of H 2 molecules, photo-electric heating due 
to external radiation (de Jong 1977) and ambipo- 
lar drift (Padoan, Zweibel & Nordlund 2000b). 



In order to derive the radiative cooling rate for 



a given species one must first know the excita- 
tion conditions in the cloud. However, the cooling 
rates also depend on the detailed structure of the 
cloud and the true net flow of energy can only 
be obtained by solving the full radiative transfer 
problem. So far the calculations have been based 
on homogeneous models in different geometries or 
spherically symmetric clouds with smooth density 
distribution. 

Because of its large abundance the CO molecules 
is usually the most important coolant of molec- 
ular gas. Goldreich & Kwan (1974) studied the 
CO cooling and made comparisons between the 
CO cooling rates and the heating associated with 
cloud collapse. The radiative transfer calculations 
were based on the Sobolev, or large velocity gradi- 
ent (LVG), approximation. Goldsmith & Langcr 
(1978) computed cooling rates due to a number 
of molecules and atomic species using the LVG 
method applied to spherical, homogeneous and 
isothermal model clouds. They covered a wide 
density range up to 10 7 cm~ 3 and a temperatures 
up to 60 K. Although 12 CO dominates the cooling 
of low density gas, the rare CO isotopes, CI and 
O2 were found to contribute a large fraction of the 
total cooling rates. At even higher densities the 
cooling rates are determined by a large number 
of species including H 2 0, various hydrides and 
molecular ions. Goldsmith & Langer considered 
several heating mechanisms and derived equilib- 
rium temperatures for typical clouds. 

Neufeld et al. (1995) continued this work by 
studying radiative cooling rates of dense molec- 
ular clouds, n = 10 3 -10 10 cm~ 3 using updated 
molecular data. The molecular gas was assumed 
to be fully shielded from external ionizing radia- 
tion. Neufeld ct al. used chemical models to com- 
pute the steady state molecular abundances that 
were used as the basis of the cooling rate calcu- 
lations. Level populations of the main cooling 
species were solved using the escape probability 
formalism (Neufeld & Kaufman 1993). Strictly 
speaking the calculations assume plane-parallel 
geometry with strong velocity gradient but the re- 
sults can be applied also to different geometries. 
That paper concentrated on the study of isother- 
mal, spherical models. 

In the present work we re-examine the cooling 
efficiency of the main cooling agents in molecular 
clouds at temperatures Tki n <100K. Compared 



with the earlier studies there are two main im- 
provements: i) The model clouds are no longer as- 
sumed to be homogeneous and, more importantly, 
the density and velocity structures are the result of 
realistic magnctohydrodynamical simulations; ii) 
The non-LTE radiative transfer problem is solved 
exactly with Monte Carlo methods. The solution 
takes fully into account the inhomogeneous den- 
sity and velocity fields of the models. 

2. The Cloud Models 

The numerical models used in this work are 
based on the results of numerical simulations 
of highly super-sonic magneto-hydrodynamic 
(MHD) turbulence, run on a 128 3 computational 
mesh, with periodic boundary conditions. 

As in our previous work, the initial density and 
magnetic fields are uniform. We apply an external 
random force, to drive the turbulence at a roughly 
constant rms Mach number of the flow. The force 
is generated in Fourier space, with power only at 
small wave numbers (1 < k < 2). The isothermal 
equation of state is used. A description of the 
numerical code used to solve the MHD equations 
may be found in Padoan & Nordlund (1999 and 
references therein). 

In order to scale the models to physical units, 
we use the following empirical Larson type rela- 
tions, as in our previous works: 

"•-"(is) • <» 

where Ai s is the rms sonic Mach number of the 
flow (the rms flow velocity divided by the sound 
speed), and a temperature T = 10 K is assumed, 
and 

{U) = 2 '° X ( 1^) ' (2) 

where the gas density n is expressed in cm~ 3 . The 
rms sonic Mach number is an input parameter 
of the numerical simulations, and can be used to 
scale them to physical units. The rms Alfvenic 
Mach number of the flow A4 a is also an input pa- 
rameter of the numerical simulations. M. a is de- 
fined as the ratio of the rms flow velocity and the 
Alfven velocity, {{B 2 ) / A-k {p))z . It determines the 
magnetic field strength, once the sonic rms Mach 
number, M s , is fixed. We refer to the turbulent 



2 



flow as super- Alfvenic when M. a > 1, while by 
cquipartition turbulence we mean M. a 1. 

In this work we use three models. They are all 
highly super-sonic, with M s ~ 10. Models B and 
C are super- Alfvenic, with M a ~ 10, while model 
A has rough equipartition of magnetic and kinetic 
energy of turbulence, with M a ~ 1- Models A 
and i? neglect the effect of self-gravity, which is 
instead included in model C. 

The physical unit of velocity in the code is the 
isothermal speed of sound, C s , and the physical 
unit of the magnetic field is C s (47r(p)) 2 (cgs). As- 
suming a kinetic temperature of Tki n =10K and a 
mean density of 320 cm~ 3 the mean field strength 
is 47.0 [iG in model A, 2.2 [iG in model B and 
2.6 /xG in model C. At Tki n =10K the rms veloc- 
ity is approximately 2.1 km/s and the linear size 
L = 6.3 pc in all three models. The turbulent ve- 
locity inside each computational cell is estimated 
as the rms velocity between neighboring cells. The 
macroscopic velocities and the thermal line widths 
are scaled according to the assumed temperature. 

In Figure 1 we show the distribution of cell den- 
sities, velocities and magnetic field strength for the 
three models. Model B has the widest range of 
densities and the largest fraction of dense cells, 
n ^10 3 cm~ 3 (top panel of Figure 1). In C the 
density distribution is skewed towards lower den- 
sities. Model A is the least inhomogeneous one, 
but the density contrast is still larger than three 
orders of magnitude. As can be seen from the mid- 
dle panel of Figure 1, differences in the velocity 
distribution between different models are insignif- 
icant. The magnetic field strengths are shown in 
the bottom panel of Figure 1 , with physical values 
obtained for the mean density n — 320 cm~ 3 and 
the kinetic temperature Tk; n = 20 K. 

3. The Calculation of the Cooling Rate 

The cooling rates are calculated individually for 
each cell in the model clouds and for several cool- 
ing species. The original MHD simulations were 
performed on a 128 3 grid but because of the high 
computational burden of the radiative transfer cal- 
culations most of these calculations were carried 
out with model clouds re-sampled into a 90 3 cell 
grid. The discretization introduces some small 
smoothing of the high density peaks but is oth- 
erwise not expected to affect the derived cooling 
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Fig. 1. — Number of cells in the three model clouds 
(A, B and C) as a function of density (a) veloc- 
ity (b) and magnetic field strength (c). The plots 
assume a mean density of 320 cm~ 3 and a temper- 
ature of T kin =20K. 

rates. 

The collisional coefficients for CO were taken 
from Flower & Launay (1985), CI-H 2 rates are 
from Schroder et al. (1991) and OI-H 2 rates from 
Jacquet et al (1992). The rate coefficients for 
collisions O2-H2 and C-2-He were provided by P. 
Bergman (1995; private communication). 
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3.1. Fractional Abundances of the Cooling 
Species 

The cooling rates are calculated for the species 
12 CO, 13 CO, C 18 0, 2 , CI and 01. According to 
Neufeld & Kaufman (1993) these are the most im- 
portant cooling species in the present parameter 
region, that is at densities below ~ 10 4 cm~ 3 and 
at temperatures T^ n <100K. Cooling rates due to 
H 2 are calculated only at Tki n =60K. 

For the CO species we assume fractional abun- 
dances [ 12 CO]/[H 2 ]=5-l(r 5 , [ 13 CO]/[H 2 ]=M(r 6 
and [C 18 0]/[H 2 ]=1-1CT 7 . These values are similar 
to those adopted by Goldsmith & Langer (1978). 
The CO abundance is lower than predicted by 
standard chemical models (Millar et al. 1997; Lee 
et al. 1997) but consistent with observation (e.g. 
Ohishi et al. 1992) which also show a significant 
variations between clouds (Harjunpaa, & Mattila 
1996). Depending e.g. on the radiation field and 
the gas temperature, chemical fractionation can 
lead to abundance variations (e.g. Warin et al. 
1996) but we shall assume constant abundances 
throughout the clouds. This is probably a good 
approximation since the inhomogeneous density 
distribution reduces differences in the radiation 
field between inner and outer parts of the cloud 
(Boissc 1990; Spaans 1996) and we do consider 
only isothermal models. 

The abundances of 2 and 01 arc not well 
known. Recent results from the ISO and SWAS 
satellites show that the oxygen abundances in 
molecular clouds have been previously overesti- 
mated and probable values are below [0 2 ]/[H 2 ]^ 
1(T 6 (Bcrgin et al. 2000; Goldsmith et al. 2000). 
On the other hand, Caux et al. (1999) found 
L1689N to be rich in OI indicating that the abun- 
dance of OI can be as high as [OI]/[H 2 ]> 10~ 5 . 
The predictions of chemical models have usu- 
ally been closer to 10~ 4 for both [OI]/[H 2 ] and 
[0 2 ]/[H 2 ] (e.g. Lee et al. 1997). Above 100 K the 
abundances are also very sensitive to the assumed 
temperature but in the temperature range consid- 
ered in the present work the abundances are es- 
sentially independent of temperature (Neufeld et 
al. 1995). We use an abundance of 10~ 5 for both 
OI and 2 . In view of the recent observational 
results the OI abundance is adequate but the 2 
abundance might be too high by more than one 
order of magnitude. However, a similar value is 



used by Goldsmith & Langer (1978) and Neufeld 
et al. (1995), which makes comparison with their 
results easier. For CI a fractional abundance value 
of 10~ 6 will be used. 

For water a fractional abundance of 1.0- 10 -6 
is assumed. Calculations are carried out sepa- 
rately for ortho- and para-water with ortho to 
para ratio 1:3. The total abundance is similar 
to the values used by Goldsmith & Langer (1978) 
and Neufeld et al. (1995). However, observations 
of quiescent gas in molecular cloud cores, Orion 
and M17SW (Snell et al. 2000a,b,c) have indi- 
cated much lower abundances [H 2 0]/[H 2 ]^ 10~ 8 . 
Similar values have been reported by Ashby et al. 
(2000). The lower values would mean that water 
is unimportant for the cooling of the clouds con- 
sidered here. Locally H 2 can be very efficient 
coolant since in outflows and hot cores its abun- 
dance can be enhanced up to [H 2 0]/[H 2 ]^ 10~ 4 
(Snell et al. 2000c) and the importance of H 2 
increases with temperature. Our models repre- 
sent, however, relatively cold and quiescent clouds. 
The H 2 cooling rates are therefore computed 
only for Tki n =60K and with relative abundance 
[H 2 O]/[H 2 ] = 1.0-10- 6 the results can be taken as 
upper limits for the actual H 2 cooling. 

Since we will examine only isothermal mod- 
els, the possible temperature dependence of the 
abundances affects only the comparison between 
models with different Tkin- At these low temper- 
atures the temperature dependence is, however, 
weak. According to models (e.g. Lee et al. 1997) 
there can be a significant dependence on the gas 
density even when photo-processes due to external 
radiation field are not considered. The effect is es- 
pecially clear for carbon. In the standard model 
presented by Lee et al. (1997) the carbon abun- 
dance increases by a factor ~10 as the gas den- 
sity decreases from 10 cm" 3 to 10 3 cm" 3 . This 
would reduce the spatial variation in the carbon 
emission. We have not included these abundance 
variations in our models. The steady state abun- 
dances are not necessarily valid for the turbulent 
medium of the MHD models where density varia- 
tions are caused by moving shock fronts and the 
predicted abundances must be treated with some 
caution. 
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3.2. Energy Levels Used in the Calcula- 
tions 

For practical reasons the number of energy lev- 
els that can be included in the calculations is lim- 
ited. The number of levels needed for an accurate 
estimate of the cooling rate depends mainly on 
the excitation as all significantly populated levels 
must be considered. The relative importance of 
the transitions is affected also by optical depth, 
as optically thick transitions contribute less to the 
total cooling rate. 

In the case of 12 CO we find that at low tem- 
peratures, Tki n < 20 K, it is sufficient to include 
levels up to J = 11. Higher transitions are more 
important at higher densities. However, at 20 K, 
the cooling rate of the J = 8 — 7 transition is ap- 
proximately two orders of magnitude below that 
of the J = 2 — 1 transition, even for densities 
n >10 3 cm~ 3 . The total contribution from higher 
transitions is insignificant. 

At higher temperatures, Tki n > 20 K, we in- 
clude the 15 lowest energy levels of 12 CO. In Fig- 
ure 2 the cooling rate in model B is plotted as 
a function of the kinetic temperature. At Tki n > 
60 K and at densities n ^10 4 cm -3 the cooling rate 
from the J = 8 — 7 transition is roughly equal to 
the rate from the J = 2 — 1 transition. However, 
the population as a function of J decreases rapidly 
and the contribution of levels close to J = 14 is 
insignificant. All models have similar ranges of 
density and the same number of levels is used also 
for the other two model clouds (A and C). 

Due to small optical depths the number of pop- 
ulated 13 CO and C 18 levels is lower than in 
12 CO. For these species we include 10 levels at 
Tkin <20K and 12 at higher temperatures. The 
number of O2 levels used in the calculations was 
the same as for CO but because of the fine struc- 
ture there are now 22 transitions between the first 
15 energy levels. For both CI and OI, 3 energy 
levels were used, with only 2 transitions between 
them. The CI lines are 3 Pi- 3 Po at 610 /im and 
3 P 2 - 3 P at 230 /im and the corresponding OI lines 
are at 145 «m and 63 «m. Cooling rates due to 
H2 O were calculated only at 60 K and only 1 1 en- 
ergy levels (20 transitions) were included. This is 
adequate due to the relatively low densities and 
temperatures. The highest included levels arc 
more than 400 K above the ground state and are 



not significantly populated. 

The inhomogeneous cloud structure of the 
MHD models increases the photon escape proba- 
bility in otherwise optically thick transitions. The 
cooling rate from lower transitions is increased 
while the population of the higher levels is de- 
creased. 
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Fig. 2. — The cooling rates, A/n, due to selected 
12 CO transitions (J =2-1, 4-3, 6-5 and 8-7) as 
a function of the kinetic temperature. The rates 
are averaged over cells in the density ranges 10 3 - 
10 4 cm -3 (upper panel) and 10-10 2 cm~ 3 (lower 
panel). The errorbars correspond to the disper- 
sions within the given density ranges. 

4. The Radiative Transfer Method 

The model cloud is divided into 90 3 cells. Each 
cell is assumed to be homogeneous and is charac- 
terized by one value of density, intrinsic linewidth 
and macroscopic velocity. Density and velocity are 
obtained directly from MHD simulations once the 
results are scaled to physical units. The intrinsic 
linewidth is the sum of thermal line broadening 
and the Doppler broadening caused by turbulent 
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motions inside the cell. The latter is estimated 
as the rms velocity difference between neighbor- 
ing cells. The components are of the same order 
but usually the turbulent line broadening is the 
larger of the two. 

The radiative transfer problem is solved with a 
Monte Carlo method (method B in Juvela 1997). 
The radiation field is simulated by a large num- 
ber of photon packages going through the cloud. 
These represent both photons entering the cloud 
from the background and photons emitted within 
the cloud. As a photon package goes through a 
cell the number of photons absorbed within that 
cell is removed from the package. At the same 
time the number of upward transitions induced 
by these photons is stored in counters. There is a 
separate counter for each cell and each simulated 
transition. After the simulation of the radiation 
field the counters are used to obtain new estimates 
of the level populations. The whole process is re- 
peated until the level populations converge (the 
relative change from an iteration to the next is less 
than 10~ 4 ). The core saturation method was used 
to speed the calculations of optically thick species 
(Hartstein & Liseau 1998; Juvela & Padoan 1999). 

The cooling rates of the cells are also calculated 
with a Monte Carlo method, using the previously 
obtained level populations. During the computa- 
tion the net flux of photons is counted for each cell 
and this is transformed into cooling rates in units 
of ergs -1 cm~ 3 . In order to study the relative im- 
portance of the transitions the net flux is counted 
for each transition separately. 

5. Results 

5.1. Cooling Rates 

The local cooling rate depends mainly on three 
parameters: the local density, kinetic temperature 
and effective column density, or effective optical 
depth. The optical depth determines the photon 
escape probability that in our case is strongly af- 
fected by the inhomogencity of the clouds. The op- 
tical depth seen by a cell varies strongly depending 
on the line of sight. This is true throughout the 
cloud, not only close to the cloud surface. The ef- 
fective optical depth depends on both the density 
and the velocity distributions since the velocity 
dispersion is always large compared with the ther- 
mal linewidth. The velocity dispersion in the mod- 




10 

N eff (H 2 ) [cm" 2 ] 



Fig. 3. — Distributions of effective hydrogen col- 
umn density, N e fj, in the three model clouds. The 
average density is 320 cm -3 and the linear size 
6.3 pc in all models. 

els is approximately 3.0 km/s, assuming a sound 
speed 0.3 km/s. 
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Fig. 4. — 12 CO cooling rates A/n as a function of 
density, n, and effective column density, N e s, in 
model B. 

In order to quantify these effects, we have cal- 
culated the effective column density, iV c ff, seen by 
each cell in the three models. We define this as 



N, 



off 



( I N(v)<f>(v)dv) (3) 



N(v) is the column density along one line of sight 
from the edge of the cloud to the cell and <j>(v) 
is the local absorption profile, both given in units 
of velocity. The integral is proportional to the 
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column density seen by the cell towards one direc- 
tion and the averaging is done over all directions. 
A^eff is calculated as a by-product of the radiative 
transfer calculations. The cooling rate depends 
on N e ff through its effect on the effective optical 
depth seen by the cells i.e. the escape probabil- 
ity. The optical depth depends in a complicated 
way on the excitation in other parts of the cloud. 
The effective column density is, on the other hand, 
a parameter of the cloud itself and describes the 
effects of the density and velocity fields indepen- 
dently of the studied molecule. For these reason 
the effective column density will be used as a sub- 
stitute for the effective optical depth. 

The distributions of N e g seen by individual cells 
are shown in Figure 3 for the three models. -/V e ff 
is smallest on the surface of the cloud and in- 
creases towards the centre depending, however, 
on the actual density and velocity fields. In mi- 
croturbulcnt case with gaussian line profiles we 
would get N c g = N/^/A-kct 2 . Here N is approx- 
imately equal to the average column density be- 
tween the cell and the cloud boundary. In our 
case we have < n >= 320cm -3 and L/2=3.1pc. 
When total velocity dispersion, some 3kms _1 at 
20 K, is used for a we obtain an effective column 
density of ~3T 20 cm" 2 km -1 s in the cloud cen- 
tre. In MHD simulations the turbulence is, how- 
ever, not random and more importantly the den- 
sity distribution is not constant. Higher values, up 
to N c ff ~ 10 21 cm~ 2 km" 1 s, are therefore reached 
in dense and velocity coherent regions. 

Figure 4 summarizes the CO cooling rate for 
model B as a function of density, n, effective col- 
umn density, N e s, and temperature. In Figure 5 
the cooling rate A(CO)/n in model B is plotted 
as a function of the effective column density N e g, 
for different density intervals and kinetic temper- 
atures Tkin- This illustrates some basic features 
seen in the case of all three model clouds. 

K/n decreases with N e g, particularly in the 
case of low temperatures and high volume den- 
sities. The behavior is the result of the dense 
cores becoming opaque. In low density regions, 
the net cooling is reduced also by the incoming flux 
from surrounding regions that, due to the higher 
density, have higher excitation temperature. At 
higher kinetic temperatures more transitions are 
contributing to the cooling rates and photon trap- 
ping cannot reduce the cooling rates to the same 
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Fig. 5. — Cooling rate A(CO)/n in model B as 
a function of the effective column density N e s, 
in clouds with kinetic temperature Tki n =10, 20 
or 60 K. Each panel shows the average cooling 
rates in three density intervals, 2.2-4.2-10 1 cm" 3 
(lowest curves), 2.2-10 2 -4.2-10 2 cm" 3 and 2.2-10 3 - 
4.2-10 3 cm" 3 (highest curves). The errorbars re- 
flect the variation of A(CO)/n in the given density 
intervals. 



extent. 

Similar effects can be seen by studying the cool- 
ing rates K/n computed separately for different 
CO transitions. Figure 6 illustrates the situation 
for model B. The increasing importance of higher 
CO transitions with increasing column density is 
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evident. In the higher density ranges the impor- 
tance of higher transitions is larger (see Figure 6). 
The increase in the cooling rate of the J = 6 — 5 
transition, for example, tends to be more rapid 
than the decrease in lower transitions. As a result, 
the total cooling rate does not drop significantly 
with N e g. 
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Fig. 7. — Cooling rate A(CO)/n from different 
CO transitions in model B, as a function of the 
effective column density N c g at three tempera- 
tures, Tkin=10K (top panel), 20 K (middle panel) 
and 60 K (bottom panel). The curves represent 
the average rate for the density interval 2.2-10 2 - 
4.2-10 2 cm" 3 . 



Fig. 6. — Cooling rate A(CO)/n from differ- 
ent CO transitions in model B, as a function of 
the effective column density iV c ff. The kinetic 
temperature is 7ki n =20K. The panels correspond 
to the density ranges: 2.2-10 1 -4.2-10 1 cm" 3 (a), 
2.2-10 2 -4.2-10 2 ckT 3 (b) and 2.2T0 3 -4.2-10 3 cm" 3 
(c). The variation in the cooling rates within the 
given density intervals is indicated by errorbars. 

Figure 7 shows similar A/n dependencies at 



three temperatures. At Tki n = 10 K the cool- 
ing rate of the transition J = 2 — 1 (and J = 
1 — 0) decreases strongly as column density ex- 
ceeds 10 20 cm -2 . The increase in the cooling rate 
from the higher transitions is unable to compen- 
sate for the loss and the total cooling rate de- 
creases as seen in Fig 5. On the other hand, al- 
ready at Tkin=20K the cooling rate by the tran- 
sition J = 4 — 3 exceeds at high column densities 
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that of the transition J = 2 — 1 and A/n levels off. 
However, the rates do not decrease significantly 
even at N cS = 5- 10 20 cnT 2 . At T kin =20K the av- 
erage optical depths averaged over the whole cloud 
are r( J = 2 — 1) w 130, r( J = 4 - 3) w 65 and 
t(J = 6 — 5) w 0.5. At 60 K the corresponding 
values are r( J = 2 — 1) w 30, r( J = 4 - 3) w 42, 
t( J = 6 — 5) w 11 and r( J = 8 - 7) w 0.4. At 
higher temperatures optical depths of individual 
transitions are lower and this reduces the effect 
that column density has on local cooling rates. 
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Fig. 8. — The cooling rate A/n in the cells of 
model B as a function of local gas density. The 
curves correspond to different values of the kinetic 
temperature, Tk; n =10, 20, and 60 K, in increasing 
order. Rates are plotted for every tenth cell in the 
model. 

5.2. Analytic Approximations 

We have computed an analytic approximation 
of the cooling rates A = A(n, T^ n , N c g). The rates 
A/n are fitted at each temperature with a function 

f(n,N eS ) = ax 
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Fig. 9. — Parameters of Eq. 4 as a function of 
the kinetic temperature of the model, Tki n , when 
fitted to the 12 CO cooling rates A(CO)/n. The 
panels correspond to the three cloud models A, B 
and C. 
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This is not the most optimal functional form for 
fitting the cooling rates but is conceptually sim- 
ple. Parameters b and d are, respectively, the 
slopes of the density and the column density de- 
pendence and the parameters c and e represent 
the non-linearity of these relationships. As a first 
approximation the dependence of the A/n on the 
density and column density is linear on the log-log 
scale. Deviations from this behavior are visible at 
low kinetic temperatures as a flattening or even 
a turnover in the density dependence. Although 
similar turnover is seen in the cooling rates of indi- 
vidual CO transitions also at higher temperatures 
(see Figure 8) the total cooling rates are monotonic 
in the present density and column density ranges. 
The turnover is simply transferred to higher den- 
sities and/or column densities. In the case of CI, 
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flattening is even more pronounced and persists to 
higher temperatures. 

For the fitting of Eq. 4 the cells of the clouds 
were divided into small density and column den- 
sity bins and the parameters of Eq. 4 were fitted 
using the average values (n), (N c s) and (A) in 
each bin. The fitting was weighted with the num- 
ber of the cells in each bin and the fit is therefore 
least reliable in the tails of the density and column 
density distributions. 

Figure 9 shows the fitted parameters as a func- 
tion of the kinetic temperature in the case of 12 CO. 
Note the decreasing dependence of cooling rate on 
column density when kinetic temperature is in- 
creased (coefficients d and e). As already stated, 
this is a natural consequence of the reduced op- 
tical depth per transition. A related effect is the 
increased importance of density at higher temper- 
atures (coefficients b and c). The first excitation 
level of CO is at ^5K and e.g. J = 5 is <~ 80 K 
above the ground state. Therefore, at T^ n ~10 K 
only a couple of levels can be populated whereas at 
2~kin ^80K a density increase can easily double the 
number of populated levels and markedly increase 
the escape probability of the emitted photons. 

The dependence on the effective column den- 
sity, N c g, is not very strong for the cooling rate 
A/n. This is true for optically thin species and 
even for 12 CO, at least for temperatures above 
10 K (see Figure 5). It is therefore reasonable 
to search for an analytical approximation for A/n 
as a function of density and kinetic temperature 
alone. The function 

log (A/n) = c lv /Iogr+(logn) 3 / 2 (c 2 + c 3 /^if 

(5) 

may be used to represent the 12 CO and the total 
cooling rate with ~20% accuracy over the studied 
density and temperature ranges, where the abso- 
lute values of A/n change by more than a factor 
of 10 4 . The parameter values obtained for the fits 
to total cooling rates are listed in Table 1. 

6. Discussion 

In the parameter range studied 12 CO is the 
main coolant. In Figure 10 we plot the cooling 
rates from model B as a function of the gas den- 
sity, for all the coolants included in our study. The 
most important difference between Tki n =20K and 
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Fig. 10. — Cooling rate, A/n, of the studied 
species as a function of the gas density in model 
B. The lower panel corresponds to kinetic temper- 
ature Tkin=20K; the upper panel to Xkin=60K. 
The rates for CI are drawn with dashed lines. 
At 60 K the sum of the computed ortho-H^O and 
para-H 2 rates is shown as a dotted line. With 
the abundances used in this paper the O2 rates 
(not shown) would be roughly equal with the CI 
rates. As the real O2 abundances are expected to 
be much lower its contribution to the cooling will 
be negligible 

60 K is the increased cooling from O at the higher 
temperature, where it provides w 10% of the to- 
tal cooling rate. Due to the lower optical depth 
per transition the 12 CO cooling rate is, at 60 K, 
considerably higher than the cooling rate of 13 CO 
and C 18 0. The sum of the computed para-H^O 
and ortho-H 2 rates are shown for Tk; n =60K. At 
this temperature the contribution from H 2 is not 
very significant. 

According to Neufeld et al. (1995) the most im- 
portant group of coolants not included in our cal- 
culations is the non-hydride molecules. Neufeld et 
al. estimate the total effect from these using the 
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formula 

L= E*„t oroi cico(0.01n,100JV) 

+ Epoiyatoimc Leo (0.008 n, 7 N) (6) 

with L(M) defined as A = Ln(H 2 )n(M) for each 
species M. According to the formula the rates are 
similar to those of CO in gas with lower volume 
density and higher column density. The average 
dipolc moment of diatomic non- hydrides (CS, NO, 
CN etc.) is approximately 10 times the dipolc mo- 
ment of the CO molecule. The radiative rates are 
therefore two orders of magnitude higher and this 
results in the first term of the equation. The sec- 
ond term for polyatomic non-hydrides was based 
on the computed cooling function of SO2 (see 
Neufeld et al. 1995). 

Using the CO cooling functions derived in this 
work, substituting N c g for N in Eq. 6 and us- 
ing the steady state fractional abundances pub- 
lished by Lee et al. (1997), 'new standard model', 
2kin=10K, n=10 3 cm -3 , we can estimate the con- 
tribution from these other molecules. With the 
CO cooling function fitted in model B we find 
that at T=10K this amounts to -10% of the CO 
cooling rates at the high density limit of our mod- 
els, n=10 4 cm -3 . The importance of these other 
molecules decreases rapidly with decreasing den- 
sity. The cooling rate decreases also with increas- 
ing kinetic temperature so that at 60 K it is less 
than 5% of the CO cooling rate. 

Although the previous estimates are very crude 
we can conclude that for our models the non- 
hydride molecules are not important, except per- 
haps in the densest and coldest cores, where they 
could provide s=s 10% of the total cooling. How- 
ever, this number is uncertain and can be altered 
significantly, for example by assuming different 



Tabic 1: Parameters of Eq. 5 obtained by fitting 
the cooling rate, A/n [10~ 25 ergs -1 ], in the three 
models. The equation 5 gives the local cooling 
rate as a function of local gas density and the ki- 
netic temperature. 
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fractional abundances. 

6.1. Comparison Between MHD Models 

Both density and velocity fields are important 
in determining the cooling rate of optically thick 
lines. We may therefore expect to see some dif- 
ferences between the MHD models although they 
are similarly inhomogeneous in both density and 
velocity space. 

Compared with the model B at Tki n =10K the 
average cooling rate, A/n, is —15% lower in model 
A. The difference increases with temperature 
and exceeds 30% at Tkin=60K. The difference 
is, however, partly due to the fact that model 
B has more high density cells with correspond- 
ingly higher cooling rates (see Figure la). In 
fact, for densities below 10 3 cm -3 the cooling ef- 
ficiency is higher in model A. At column density 
N c g wlO 20 cm -2 and volume density 10 2 cm -3 the 
cooling rates (A/n) in model A are a few of percent 
larger than in model B. At somewhat higher den- 
sities, n w 3-10 3 cm -3 , the rates are below those 
of the model B (at 10 K by some 10% but at 60 K 
by only -1%). 

The net cooling rate (A/n) in model C is 4% 
smaller than in model B. As a function of density 
and column density the rates are between those of 
the models A and B and usually closer to the rates 
in model A. The model C included self-gravity 
which should affect the structure of the high den- 
sity concentrations. This does not, however, show 
clearly in the overall density distribution (see Fig- 
ure 1) and even at n ~10 3 cm -3 the cooling rates 
are quite similar to those found in the other mod- 
els. Same trend continues up to highest densities, 
n <~10 4 cm -3 , with model B having cooling rates 
5-10% in excess of the other two models. 

The cooling function should be sensitive to the 
cloud structure. In the range n =10 1 -10 3 cm -3 
and N =10 19 -10 21 cm -2 the differences in the 
12 CO cooling efficiency in the three MHD models 
are, however, below 20% and from this point of 
view the models do not differ radically from each 
other. All have similar types of inhomogeneous 
density structures and the overall velocity disper- 
sions are also similar. Furthermore, the dumpi- 
ness increases the photon escape probability and 
reduces any effects caused by large optical depths. 
The differences between MHD models would be 
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enhanced e.g. in clouds with larger size, higher 
density or lower turbulence. 

6.2. Average cooling rates 

In the previous chapters we have studied the lo- 
cal cooling rate as a function of the local gas prop- 
erties, density n and kinetic temperature T kin , and 
the general environment described by the effective 
column density N e g. Since earlier theoretical stud- 
ies and observers tend to concentrate on the av- 
erage properties of the clouds we will now discuss 
the average cooling rate which is also proportional 
to the net cooling of the cloud as a whole. 
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Fig. 11. — The cooling rates averaged over model 
B, as a function of the kinetic temperature of the 
model. The mean density is < n >=320cm~ 3 . 

Figure 11 shows the cooling rates in model B 
averaged over the model volume. Note the steep 
temperature dependence of the 01 emission. The 
actual ratios between different species depend on 
the assumed abundances and, as already men- 
tioned, for species other than 12 CO the cooling ef- 
ficiency scales almost linearly with the abundance. 

We now study more closely the 12 CO cooling. 
The cooling rate depends strongly on the density 
(Figure 8). In an inhomogeneous cloud this means 
that most of the cooling power comes from regions 
denser than the average. The average cooling rate 
is also higher than the local rate at the density 
equal to the mean density of the cloud. This is true 
as long as the dense cores are not optically very 



thick. When the cores become opaque to most 
of the cooling lines the cooling efficiency will not 
increase any more with increasing density and the 
total cooling efficiency will be reduced. 
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Fig. 12.— The 12 C0 cooling in model B at three 
kinetic temperatures. The longer curves show 
the local cooling rates as a function of the lo- 
cal gas density. The mean density of the model 
was 320 cm~ 3 . The filled squares show the cool- 
ing rates averaged over the model, as a function of 
the mean density of the model. The models were 
obtained by scaling the mean density of model B 
to 160, 320, 640, 1280 and 2560 cm" 3 

Figure 12 shows both the local cooling rate and 
the average cooling rate of 12 CO, as a function 
of density for three different kinetic temperatures. 
The local rates are obtained from the individual 
cells of the model B as the cooling rate of a cell 
divided by the gas density in the cell. The entire 
model cloud provides the point at 320 cm -3 for the 
curve showing the average cooling rate < A > / < 
n >. Here the cooling rate and the density are 
averaged separately over the volume of the model. 
The average cooling rate is proportional to the net 
cooling rate of the entire model (ergs^ 1 ). Other 
points for the average cooling rate are from models 
which were obtained by scaling the mean density 
of model B to 160, 640, 1280 and 2560cm~ 3 . 

The rates follow the expected behavior outlined 
above. At lower densities the average rates exceed 
the local rates, with a ratio approaching a factor 
of ten. Even in the < n >=320cm~ 3 model the 
CO rate A/n levels off at the highest densities. 
In models with high mean density this leads to 
a turnover and eventually the average rate drops 
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below the curve for the local cooling rate. This 
means that gas below the mean density is getting 
more and more important for the net cooling of 
the cloud. For 10 K models the curves meet be- 
low 10 3 cm~ 3 . At higher temperatures the opti- 
cal depth effects are again reduced and at 60 K 
the average rate remains above the local rate until 
close to 10 4 cm -3 . This is perhaps the most strik- 
ing result of the present work. Molecular clouds 
described by turbulent density and velocity field 
as modeled here, can radiate 10 times more effi- 
ciently than a uniform cloud with density equal to 
the mean density of our model. As a consequence, 
their thermal balance requires a ten times larger 
heating source. 

In order to study separately the effects of the 
density and the velocity inhomogeneities we have 
considered two further models that are based on 
model B. The constant density model, CD, has 
constant density 320 cm~ 3 (equal to the mean den- 
sity in the model B) but the velocity structure of 
the model B. The constant velocity model, CV, 
has the same density structure as the model B but 
no macroscopic velocity field. Only 12 CO cooling 
rates were calculated since, due to high optical 
depth, these are most likely to show differences 
between the models. 

Compared with model B the effective column 
densities are higher in CV and CD, since the in- 
homogeneity of either densities or velocities is re- 
moved. This gives the appearance of increased 
A/n when plotted against N e g. This is, however, 
only due to the shift of the N e s axis and the total 
cooling rate of the clouds is reduced. This is natu- 
ral since the escape probability of emitted photons 
is smaller. For model CV the total 12 CO cool- 
ing rate is reduced by ~10% at Tk; n =10K. The 
effect reduces at higher temperatures as the op- 
tical depths of individual transitions are reduced. 
For model CD the drop is ^50% at 10 K and in- 
creases with temperature. The difference between 
the original model and CD is not due to radia- 
tive transfer effects but is rather a direct conse- 
quence of the steep density dependence of the cool- 
ing function. In an inhomogeneous cloud most of 
the cooling is provided by regions with density well 
above the average value. For optically thin emis- 
sion the ratio of average cooling rate and cooling 
rate at the mean density can be derived directly 
based on the cell density distribution (Figure 1) 



and the density dependence A ~ n 2 . The ratios 
are 3.9, 8.0 and 8.9 for the models A, B and C, 
respectively. 

Since density peaks are very important for the 
cooling we must make sure that the results are 
not affected by the limited spatial resolution. The 
CO cooling rates were compared in three variants 
of the model B where the cloud was divided into 
cither 128 3 , 90 3 or 48 3 cells. The plots of local 
cooling rates against local density were almost in- 
distinguishable although at lower resolutions the 
total range of densities was reduced. A change in 
the discretization causes changes in the cell densi- 
ties. In the plot this means only a small displace- 
ment along the A(n) curve which itself remains 
unchanged. As the cell size is increased some of 
the velocity dispersion between cells ('macrotur- 
bulence') is transformed into turbulent velocity in- 
side the cells ('microturbulcnce') but this did not 
produce any noticeable effects on the cooling rates. 
At Tki n =20K the total cooling rates of the mod- 
els were within ^3% of each other. This shows 
that the selected resolution was sufficient to cap- 
ture the effects of density and velocity inhomo- 
geneities. Our models represent quiescent diffuse 
clouds. In models with more developed cloud cores 
and steeper density gradients the spatial resolu- 
tion should be correspondingly higher. 

6.3. Comparison with Earlier Work 

Most previous calculations of molecular line 
cooling rates have focused on specific objects (for 
example dense cores) and are therefore not useful 
for comparison with the present work. In particu- 
lar, little has so far been published on the subject 
in connection with inhomogeneous clouds. In the 
following we look at some of the differences be- 
tween our work and the results obtained by Gold- 
smith & Langer (1978) and Neufeld et al. (1995). 
Their results apply to clouds with continuous and 
smooth density distributions. Differences will be 
caused by two factors. Firstly, our MHD mod- 
els contain a range of densities and due to steep 
density dependence (e.g. Figure 8) total radi- 
ated energy will be larger than in a homogeneous 
cloud with equal mean density (see previous chap- 
ter). Secondly, excitation and escape probability 
of emitted photons will be affected by the density 
and velocity inhomogeneities. 

In Figure 13 we plot the Goldsmith & Langer 
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Fig. 13. — Comparison of our 12 CO and CI cool- 
ing rates with Goldsmith & Langer (1978). The 
solid curves show our local cooling rates in model 
B as a function of the local gas density. The error- 
bars indicate the variation within the density bins. 
The solid squares show our estimates for the av- 
erage cooling as a function of the average density 
of the model cloud (see text). The dashed curves 
are 12 CO and CI rates taken from Goldsmith & 
Langer (1978) 

(1978) cooling rates for 12 CO and CI together with 
our average rates at a few densities (solid squares). 
In the same figure the local cooling rates for the 
n = 320 cm~ 3 model are also shown (solid curve). 
This should be compared only with those points 
on the other curves that correspond to the same 
density. 

We will first compare the CO rates. For the 
LVG model used by Goldsmith & Langer the rel- 
evant parameter is the CO abundance divided by 
the velocity gradient, Xqo/ (dV/dr). In our mod- 
els a corresponding parameter can be calculated 
using the average velocity dispersion of the cloud 
and the linear size of the model, L. Choosing 
models where these parameters agree ensures that 



the average properties of the LVG and the MHD 
models will be similar. Differences will be caused 
by differences at smaller scales i.e. mostly by the 
density and velocity inhomogeneity of our models. 
In the MHD simulations the velocity and density 
fields arc not completely random and depending 
on the line of sight the average density and column 
density can differ significantly from the values av- 
eraged over the entire volume of the cloud. 

In the following we will use the FWHM of 
the one-dimensional velocity distribution as the 
measure of the velocity dispersion in our mod- 
els. The value is calculated as the average over 
the whole cloud. In model B the FWHM of the 
velocity dispersion is 2.90 km s -1 at Tk; n =10K. 
With abundance 5-10~ 5 and cloud size 6.25 pc 
we get a column density per velocity unit of 
1.08T0 -4 pc km _1 s. In Figure 13 the model 
is compared with Goldsmith & Langer calcula- 
tions for X/(dT//dr)=10- 4 pckm- 1 s. At 40K 
the velocity dispersion in the MHD model was 
higher, FWHM=5.77kms -1 , corresponding to 
the higher speed of sound. The column density 
is 5.4-10~ 5 pckm _1 s and the corresponding value 
of the Goldsmith & Langer model shown in Fig- 
ure 13 is 4-10 -5 pckm _1 s. Goldsmith & Langer 
note that their results are not sensitive to the as- 
sumed velocity gradient and in our case similar 
conclusion can be drawn from the weak column 
density dependence in Figure 5. Therefore, an 
absolute equality of the column density parame- 
ters is not crucial and the qualitative result of the 
comparison would remain the same even if our 
column density values were scaled by a factor of 
two. 

At 10 K our local 12 CO cooling rate at n =320 cm -3 
agrees with the predictions of Goldsmith & Langer 
(see Figure 13). The point is still on the linear 
portion of the curve i.e. radiative trapping is not 
yet significant. The average cooling rate is, how- 
ever, two times higher than either the Goldsmith 
& Langer value or our local cooling rate at that 
density. This is due to the fact that the density 
dependence of the cooling rate A is steeper than 
n 1 ' . In the case of a wider distribution of den- 
sity values the average cooling rate increases even 
when the mean density remains unaltered. 

At 40 K the effects of optical depth are greatly 
reduced and in the model with mean density 
320 cm -3 the local cooling rate is almost a linear 
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function of local density. This increases the differ- 
ence to Goldsmith & Langer results which are at 
n = 10 4 cm -3 a factor of six below our rate. At 
10 K the average cooling rate A/n was a decreas- 
ing function of density but at 40 K it increases up 
to < n >^ 10 3 cm -3 . In the model with mean 
density < n >=320cm -3 the local cooling rate 
exceeds the volume averaged rate by a factor of 
four. 

In Figure 13 the errorbars indicate the varia- 
tion in local cooling rate, A/n, within each den- 
sity bin (see also Figure 8). Although this does 
include some noise from the Monte Carlo calcu- 
lations most of the variation is due to radiative 
transfer effects. The cells are in different environ- 
ments (e.g. cloud centre vs. cloud surface) and 
this affects the cooling. At 40 K the local cooling 
rate is rather uniform while at 10 K there is wider 
scatter, especially at high densities. This shows 
again how an increasing kinetic temperature re- 
duces the effects of optical depth. 

For CI the abundance used in this paper was 
10 -6 i.e. a factor of ten lower than the value used 
by Goldsmith & Langer. The relevant parame- 
ter is, however, again Xqi/ (dV/dr). None of the 
results published by Goldsmith & Langer corre- 
spond exactly to the parameters of our models and 
therefore we have rescalcd our abundance value 
so that our model can be compared directly with 
their Figure 7. Another difference is caused by the 
collisional coefficients. Goldsmith & Langer used 
for CI-H2 collisions the Launay & Roueff (1977) 
CI-H rates divided by ten. In this paper we have 
used rates given by Schroder et al. (1991) and 
these result in a significantly higher CI cooling. 
However, in Figure 13 we have derived the CI cool- 
ing using the scaled Launay & Roueff (1977) co- 
efficients. After these modifications the predicted 
cooling rates agree at low column densities and the 
mean density and the average column density per 
velocity interval are identical to the values in the 
Goldsmith & Langer model. 

The main features are similar as in the case of 
CO. Our local CI rate A/n increases up to den- 
sity n ~ 10 4 cm -3 where the Goldsmith & Langer 
curves are already clearly decreasing. At 40 K the 
turnover is shifted further to a higher density. The 
average CI rates at five mean densities between 
160 cm -3 and 2560 cm -3 are shown in the same 
figure. In the model < n >=320cm -3 the average 



rate exceeds the local cooling rate at this density 
by a factor of ~3. The ratio is the same at both 
10 K and 40 K. When the mean density exceeds 
1000 cm -3 the average rates drop close to Gold- 
smith & Langer predictions. 

Compared with CO the main difference is that 
the kinetic temperature has very little effect on 
the shape of the CI curves. At 40 K the CO curve 
has become almost linear while the CI has simi- 
lar flattening as at 10 K. The effect is visible also 
in the average rates. While for CO A/n changes 
from a decreasing function to an increasing one, 
little change is seen in the CI curves. There are 
only three populated CI levels and the effect ki- 
netic temperature can have on the optical depth 
of the transitions is correspondingly smaller. The 
second excitation level of CI is more than 60 K 
above the ground state while CO has already five 
energy levels below this. 

In comparison with Neufeld et al. (1995) some 
of the differences are due to the abundances. For 
example, the chemical models of Neufeld et al. 
(1995) predict OI and O2 fractional abundances 
slightly below 10 4 . The values adopted in this pa- 
per are lower by almost a factor of ten, but OI 
still provides a few per cent of the total cooling 
rate at the highest density and temperature in our 
models. The optical depth of most species is suf- 
ficiently low, so that the cooling rate depends lin- 
early on the abundance. For example, reducing 
the abundance of O2 by a factor of ten decreases 
the cooling rate by the same factor, the accuracy 
of the linear relation being «3% at 10 K and better 
than 0.5% at 60 K. On the other hand, the 12 CO 
cooling rate decreases with N e g, especially at high 
column densities. Therefore, a change in the 12 CO 
abundance has a relatively small effect on the cool- 
ing efficiency. We checked this by reducing the 
12 CO abundance by 50% in the model A. Assum- 
ing T kin =10K the value of A/n is reduced by 47% 
for column densities N e g ~ 10 19 cm -2 roughly in- 
dependent of volume density. However, when the 
column density reaches 5 • 10 20 cm -2 the reduction 
in A/n is no more than 20%. 

Neufeld et al. (1995) give cooling rate as a func- 
tion of column density measure, N — Gn/\dv/dx\, 
i.e. number density divided by the velocity gradi- 
ent or, in the absence of a gradient, by the velocity 
dispersion divided by distance. This is essentially 
the number of hydrogen atoms per velocity unit, 
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Fig. 14. — Cooling rate of 12 CO averaged over 
the volume of the model B. In the upper panel 
the filled squares joined by lines are our results 
for cloud mean densities < n >=160cm (lowest 
curve), 640cm -3 and 2560cm -3 (thick line). The 
corresponding rates computed from the formulae 
presented by Ncufeld et al. (1995) are shown with 
solid lines and dotted lines correspond to column 
densities two times or half the estimated values of 
N. In the lower panel the corresponding curves are 
plotted as the function of density for temperatures 
10, 20 and 80 K. The dashed lines are Neufeld et 
al. predictions for density 3< n >. 



in units of cm -2 km -1 s. There is, however, an 
additional factor, G, which depends on the geom- 
etry. For plane parallel medium with large velocity 
gradient the value is G = 1 and in the centre of a 
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Fig. 15. — Comparison of the cooling rates A com- 
puted from the analytical representation given by 
Ncufeld et al. (1995; dotted lines) and the rates in 
our model A (solid curve) . The curves correspond 
to density n ~ 10 3 cm -3 . The solid curve is com- 
puted from the analytical fit of Equation 4 and is 
drawn only for the actual column density interval 
present in the model. 



static sphere G ~ 0.5. The effect of the abundance 
is included in the parameter N and its value need 
not be known separately. 

In the upper panel of Figure 14 we plot the aver- 
age 12 CO rate for model _B, with the mean density 
scaled to 160, 640, and 2560cm -3 . The parame- 
ter N was estimated as GnL/a with G = 1 and a 
equal to the rms velocity dispersion in the MHD 
model at given temperature. The lower three lines 
show predictions by Neufeld et al. (1995) and the 
dotted lines show the rates with N multiplied by 
two or by one half. The lower panel shows the 
corresponding rates as function of the mean cloud 
density. The main difference is again due to the 
density inhomogeneity of MHD models. Gas with 
density above the average value provides most of 
the cooling power and the cooling rate is compara- 
ble to that of a much denser homogeneous cloud. 
In the low density models the Neufeld et al. for- 
mula would give a fairly good prediction of the 
cooling rate if density value n <~ 3 < n > were 
used. On the other hand, rate could be raised 
only slightly by lowering N. 

Next we will discuss again the local cooling 
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rates. In Figure 15 our 12 CO cooling efficiencies, 
A, are shown together with the results of Neufeld 
et al. (1995, Table 3) as a function of N. Our 
curves correspond to the analytic approximation 
of the local cooling rate in model B (see Sect. 5.2). 
For the plot we must determine a relation between 
our parameter N e g, which is a quantity integrated 
over the local absorption profile, and N used by 
Neufeld et al. For gaussian lines we have relation 
iV e fF «0.3 N/a (see Section 5.1) while for plane- 
parallel flow (geometrical factor G = 1) N = N/a. 
This gives a conversion N c g ~ 0.37V and this re- 
lation is used to plot our results on the N scale 
in Figure 15. The scaling between N and N e g is 
only approximate and depends on the LVG model 
assumed. The curves may therefore be shifted in 
the horizontal direction. However, the shift should 
correspond to no more than a factor of two change 
in N. 

The most important difference is the marked 
flattening of the column density dependence that 
is seen in our models at higher kinetic tempera- 
tures. This is the same effect as seen e.g. in Fig- 
ures 12 and 13 and we interpret this as the conse- 
quence of varying excitation conditions that allow 
efficient cooling even in the dense gas. Low den- 
sity gas radiates in low transitions while cooling 
in dense cores takes place mainly through higher 
transitions. At higher kinetic temperatures more 
transitions can contribute to the cooling and there 
is a clear difference between subthermally excited 
low density gas and cores that are close to thermal- 
ization. Different parts of the cloud are therefore 
decoupled not only due to velocity differences but 
also because they radiate mostly in different tran- 
sitions. The difference between our results and 
those in Neufeld et al. (1995) decreases with de- 
creasing temperature, i.e., at low kinetic temper- 
atures our models behave more like homogeneous 
clouds. At low temperatures the number of excited 
levels is small and the remaining transitions have 
higher optical depth. In Figure 15 we show curves 
only for density n ~10 3 cm~ 3 but the results are 
qualitatively similar even at other densities. 

The detected differences between inhomoge- 
neous and homogeneous cloud models can be at- 
tributed to two effects. Firstly, density and ve- 
locity inhomogeneities increase the photon escape 
probability. Both the excitation and the escape 
probability of photons are determined by the op- 



tical depths towards different lines of sight. In the 
center of a homogeneous cloud the optical depth 
is the same towards all directions. An inhomoge- 
neous cloud with the same average optical depth 
has always a higher escape probability, since the 
escape probability is proportional to the average 
of exp(-r) and not to the average of r. This 
leads to higher escape probability but generally 
also to lower excitation. Secondly, in the consid- 
ered models A/n is still an increasing function of 
density and therefore density variations tend to in- 
crease the total cooling power. The slope of A/n 
is, of course, determined by the radiative transfer 
and the two effects are closely interrelated. The 
density distribution was seen to affect the photon 
escape probability also indirectly through excita- 
tion. Some excitation levels are populated only 
in the densest regions and this leads to a partial 
decoupling between the dense cores and the sur- 
rounding low density gas. 

7. Summary and conclusions 

We have studied the radiative cooling of molec- 
ular gas at temperatures Tki n =10-80K and densi- 
ties n < 10 4 cm~ 3 , based on three-dimensional 
MHD calculations of the density and velocity 
structure of interstellar clouds. The models have 
been scaled to linear sizes ~6pc and mean densi- 
ties in the range 160-2560 cm -3 . The cooling rates 
for isothermal clouds were computed by solving 
the radiative transfer problem with Monte Carlo 
methods. 

We find that: 

• Inhomogeneous density and velocity fields 
reduce photon trapping and thus increase 
the cooling rates. In comparison with homo- 
geneous cloud models the MHD models are 
much less affected by optical depths effects. 
This is especially true at kinetic tempera- 
tures T ki „ >60K. 

• There is a clear difference between the den- 
sity dependence of local cooling rates and the 
density dependence of cooling rates averaged 
over entire clouds. 

• At low to intermediate densities most of 
the cooling power is provided by clumps 
with densities above the average gas den- 
sity. The average cooling rate for a model 
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with given mean density can be as much as 
an order of magnitude larger than the local 
cooling rate at the same density. In mod- 
els with higher average density and lower 
temperature the differences arc smaller since 
the densest parts of clouds become optically 
thick. 

• Compared with earlier models (Goldsmith & 
Langer (1978); Neufeld ct al. (1995)) our lo- 
cal cooling rates differ mainly at large den- 
sities where our rates are higher due to re- 
duced photon trapping. The volume aver- 
aged cooling rates are higher than in the ear- 
lier models typically by a factor of few. At 
higher temperatures (T > 40 K) the differ- 
ence can be almost one order of magnitude. 

• For the MHD models, the absence of the 
macroscopic velocity field would reduce the 
cooling by up to 10%. 

• The absence of density fluctuations would 
reduce cooling by ~50% at 10 K. This is 
caused mainly by the density dependence of 
the cooling rates and the radiative transfer 
effects are less important. At high tempera- 
tures (>80 K) the difference to homogeneous 
models approaches a factor of ten. 

• 12 CO is clearly the most important coolant 
over the whole parameter range studied. 

• At low temperatures 13 CO is the the sec- 
ond most important coolant, after 12 CO. At 
temperatures Tki n >60K it is exceeded by 
01, which can provide more than 10% of the 
total cooling (assuming a relative abundance 
~10~ 5 ). 

• In view of the recently observed very low O2 
and H2O abundances these species are unim- 
portant for the cooling of the type of clouds 
studied in this paper. 
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